﻿using System;
using System.Collections.Generic;
using System.Text;

namespace dllLink
{
    // Classe des estimateurs des paramètres d'un
    // modèle de Vasicek par maximum de log-vraisemblance
    internal class MaxLogLikEstimator : VasicekEstimator
    {
        private string mvName;
        public MaxLogLikEstimator(string theName, DateTime thePricingDate, TimeSpan theHistoryDepth)
        {
            mvName = theName;
            mvPricingDate = thePricingDate;
            mvHistoryDepth = theHistoryDepth;
        }


        // retourne un tableau de taille 3 contenant les estimations des
        // 3 paramètres du modèle de Vasicek (a, b et sigma : drt = a(b-rt)dt +sigma*dWt)
        public override double[] mComputeParams()
        {
            double[] myRes = new double[3];
            double[] myData;

            double myRx=0.0;
            double myRxx = 0.0;
            double myRy = 0.0;
            double myRyy = 0.0;
            double myRxy = 0.0;
            double myB = 0.0;

            double myLag = Model.mvDt; //pas de temps entre deux données, 1 jours
            
            double myalpha, myEstimSigma2, myNNu2;

            // on récupère les données de taux
            myData = Cotation(mvName, mvPricingDate, mvHistoryDepth);
            int myDepth = myData.Length;
            myRes[1]=0.0;
            for(int i=1;i<myDepth;i++)
            {
                myRx += myData[i-1];
                myRxx += myData[i-1]*myData[i-1];

                myRy += myData[i];
                myRyy += myData[i]*myData[i];

                myRxy += myData[i - 1] * myData[i];
                myB += myData[i - 1]/myDepth; 
            }

            
            myB += myData[myDepth - 1]/myDepth;
            myDepth -= 1;

            //b sous P
            myRes[1] = (myRy * myRxx - myRx * myRxy) / (myDepth * (myRxx - myRxy) - (Math.Pow(myRx, 2.0) - myRx * myRy));
            
            //a
            myNNu2 = myDepth * Math.Pow(myRes[1], 2.0);
            myRes[0] = -Math.Log((myRxy - myRes[1] * myRx - myRes[1] * myRy + myNNu2)
                    / (myRxx - 2 * myRes[1] * myRx + myNNu2)) / myLag;
            //sigma
            myalpha = Math.Exp(-myRes[0] * myLag);
            
            myEstimSigma2 = (myRyy - 2 * myalpha * myRxy + Math.Pow(myalpha, 2.0) * myRxx 
                - 2 * myRes[1] * (1 - myalpha) * (myRy - myalpha * myRx)
                + myNNu2 * Math.Pow((1 - myalpha), 2.0)) / myDepth;

            myRes[2] = Math.Sqrt(myEstimSigma2 * 2 * myRes[0] / (1 - Math.Pow(myalpha, 2.0)));

            //b sous Q
            myRes[1] = myB;


            return myRes;
        }

    }
}

